Trabajo 1 TAE - Ingresos
Equipo
13/2/2020
Librerías:
library(tidyverse)
library(data.table)
library(readxl)
library(GGally)
library(corrplot)
library(leaflet)
library(raster)
library(mltools)
library(highcharter)Lectura de datos:
raw_encuesta <- fread("databases/calidad_vida_ok.csv", encoding = "UTF-8") %>%
as_tibble()
# Unificación de las estadísticas de los barrios "ALTAVISTA" y "CABECERA ALTAVISTA"
raw_encuesta <- raw_encuesta %>% mutate(
encuesta_calidad.barrio = case_when(
encuesta_calidad.barrio == "CABECERA ALTAVISTA" ~ "ALTAVISTA",
encuesta_calidad.barrio == "CIUDADELA NUEVO OCCIDENTE" ~ "CABECERA URBANA CORREGIMIENTO SAN CRISTÓBAL",
encuesta_calidad.barrio == "AREA DE EXPANCION SAN CRISTOBAL" ~ "ÁREA DE EXPANSIÓN SAN CRISTÓBAL",
encuesta_calidad.barrio == "CABECERA SAN CRISTÓBAL" ~ "CABECERA URBANA CORREGIMIENTO SAN CRISTÓBAL",
encuesta_calidad.barrio == "SAN CRISTOBAL" ~ "CABECERA URBANA CORREGIMIENTO SAN CRISTÓBAL",
encuesta_calidad.barrio == "PROGRESO Nº 2" ~ "PROGRESO",
TRUE ~ encuesta_calidad.barrio
)
)
diccionario <- read_excel("databases/diccionario_ecv_2016.xlsx")Número de encuestados:
raw_encuesta %>%
count(`encuesta_calidad.barrio`, sort = T) %>%
arrange(n)## # A tibble: 296 x 2
## encuesta_calidad.barrio n
## <chr> <int>
## 1 EL CARMELO 4
## 2 GUAYAQUIL 4
## 3 LA ILUSIÓN 6
## 4 SAN JOSÉ 10
## 5 BOQUERÓN 15
## 6 EL PATIO 19
## 7 BARRIO COLOMBIA 21
## 8 DESCONOCIDO 21
## 9 YARUMALITO 22
## 10 EL UVITO 25
## # ... with 286 more rowsCreación de la base de datos:
encuesta_ingresos<-raw_encuesta%>%
mutate(proporcion_ingreso_extra=case_when(encuesta_calidad.p_71 == 1 ~1,
TRUE~ 0))%>%
mutate(proporcion_pension=case_when(encuesta_calidad.p_71 != 0~1,
TRUE~ 0))%>%
mutate(subsidio=case_when(encuesta_calidad.p_93 != 0~1,
TRUE~ 0))%>%
mutate(salario=case_when(encuesta_calidad.p_121 != 0~1,
TRUE~ 0))%>%
mutate(arriendo=case_when(encuesta_calidad.p_227 != 0~1,
TRUE~ 0))%>%
mutate(total_gastos=case_when(encuesta_calidad.p_232 != 0~1,
TRUE~ 0))%>%
mutate(gasto_servicios=case_when(encuesta_calidad.p_240 != 0~1,
TRUE~ 0))%>%
mutate(gasto_hogar=case_when(encuesta_calidad.p_254 != 0~1,
TRUE~ 0))%>%
mutate(encuesta_calidad.barrio = str_replace(encuesta_calidad.barrio, "ANDALUCIA", "ANDALUCÍA") %>%
str_replace("Nº 2", "NO.2") %>%
str_replace("Nº 1", "NO.1") %>%
str_replace("Nº 3", "NO.3") %>%
str_replace("AREA EXPANSION", "ÁREA DE EXPANSIÓN") %>%
str_replace("EXPANCION", "EXPANSIÓN") %>%
str_replace("AREA", "ÁREA") %>%
str_replace("BOMBONA", "BOMBONÁ") %>%
str_replace("LA ASOMADERA", "ASOMADERA") %>%
str_replace("BELALCAZAR", "BELALCÁZAR") %>%
str_replace("CALAZANS", "CALASANZ") %>%
str_replace("COLON", "COLÓN") %>%
str_replace("MIRA FLORES", "MIRAFLORES") %>%
str_replace("BARRIO FACULTAD DE MINAS", "FACULTAD DE MINAS") %>%
str_replace("CABECERA SAN ANT DE PR.", "SAN ANTONIO DE PRADO") %>%
str_replace("CARLOS E RESTREPO", "CARLOS E. RESTREPO") %>%
str_replace("URQUITA", "URQUITÁ") %>%
str_replace("LOS CERROS EL VERJEL", "LOS CERROS EL VERGEL") %>%
str_replace("CAYCEDO", "CAICEDO") %>%
str_replace("VALDES", "VALDÉS") %>%
str_replace("CERRO EL VOLADOR", "B. CERRO EL VOLADOR") %>%
str_replace("MOSCU", "MOSCÚ") %>%
str_replace("JOSELA", "JOSÉ LA") %>%
str_replace("JOSE", "JOSÉ") %>%
str_replace("EL YOLOMBO", "YOLOMBO") %>%
str_replace("PIEDRAS BLANCAS", "PIEDRAS BLANCAS - MATASANO") %>%
str_replace("BASILIA", "BRASILIA") %>%
str_replace("VILLA TINA", "VILLATINA") %>%
str_replace("LILIAM", "LILLIAM") %>%
str_replace("BOLIVAR", "BOLÍVAR") %>%
str_replace("CORREGIMIENTO PALMITAS", "PALMITAS SECTOR CENTRAL") %>%
str_replace("INES", "INÉS") %>%
str_replace("FE", "FÉ") %>%
str_replace("LUCIA", "LUCÍA") %>%
str_replace("SABIO", "SAVIO") %>%
str_replace("BERMEJAL- LOS ÁLAMOS", "BERMEJAL-LOS ÁLAMOS") %>%
str_replace("BOLÍVARIANA", "BOLIVARIANA") %>%
str_replace("EL NOGAL - LOS ALMENDROS", "EL NOGAL-LOS ALMENDROS") %>%
str_replace("JUAN XXIII - LA QUIEBRA", "JUAN XXIII LA QUIEBRA") %>%
str_replace("PROGRESO Nº 2", "EL PROGRESO") %>%
str_replace("MARIA", "MARÍA") %>%
str_replace("PLAYÓN", "PLAYON") %>%
str_replace("EL SOCORRO / LA GABRIELA", "EL SOCORRO") %>%
str_replace("FÉRRINI", "FERRINI") %>%
str_replace("LA CANDE LARIA", "LA CANDELARIA") %>%
str_replace("EL PLAYON", "PLAYÓN") %>%
str_replace("IGUANA", "IGUANÁ") %>%
str_replace("MARÍA CANO - CARAMBOLAS", "MARÍA CANO-CARAMBOLAS") %>%
str_replace("DE ABURRA", "DEL ABURRÁ") %>%
str_replace("ALTAVISTA CENTRAL", "ALTAVISTA SECTOR CENTRAL") %>%
str_replace("SECTOR CENTRAL", "CENTRO ADMINISTRATIVO") %>%
str_replace("ALTAVISTA CENTRO ADMINISTRATIVO", "ALTAVISTA SECTOR CENTRAL") %>%
str_replace("SANTA ELENA CENTRO ADMINISTRATIVO", "SANTA ELENA SECTOR CENTRAL") %>%
str_replace("PALMITAS CENTRO ADMINISTRATIVO", "PALMITAS SECTOR CENTRAL") %>%
str_replace("PROGRESO", "EL PROGRESO")
)
db_ingresos<-encuesta_ingresos%>%
group_by(encuesta_calidad.barrio, encuesta_calidad.comuna)%>%
summarize(n=n(),
proporcion_ingreso_extra=sum(proporcion_ingreso_extra==1,na.rm=TRUE)/sum(!is.na(proporcion_ingreso_extra),na.rm=TRUE),
proporcion_pension=sum(proporcion_pension==1,na.rm = TRUE)/sum(!is.na(proporcion_pension),na.rm=TRUE),
subsidio=sum(`encuesta_calidad.p_93`, na.rm=TRUE)/sum(subsidio==1,na.rm = TRUE),
salario=sum(`encuesta_calidad.p_121`, na.rm=TRUE)/sum(salario==1,na.rm = TRUE),
arriendo=sum(`encuesta_calidad.p_227`,na.rm=TRUE)/sum(arriendo==1,na.rm = TRUE),
total_gastos=sum(`encuesta_calidad.p_232`,na.rm=TRUE)/sum(total_gastos==1,na.rm = TRUE),
gasto_servicios=sum(`encuesta_calidad.p_240`,na.rm=TRUE)/sum(gasto_servicios==1,na.rm = TRUE)
)%>%ungroup()%>%
mutate(subsidio= replace_na(subsidio, 0))%>%
mutate(salario= replace_na(salario, 0))%>%
mutate(arriendo= replace_na(arriendo, 0))%>%
mutate(total_gastos= replace_na(total_gastos, 0))%>%
mutate(gasto_servicios= replace_na(gasto_servicios, 0))
# Eliminamos el barrio Desconocido
db_ingresos <- db_ingresos %>%
filter(encuesta_calidad.barrio != "DESCONOCIDO")
db_ingresos[complete.cases(db_ingresos),]## # A tibble: 308 x 10
## encuesta_calida~ encuesta_calida~ n proporcion_ingr~
## <chr> <chr> <int> <dbl>
## 1 AGUAS FRÍAS ALTAVISTA 126 0.00794
## 2 ALDEA PABLO VI POPULAR 426 0.00235
## 3 ALEJANDRÍA EL POBLADO 709 0.00987
## 4 ALEJANDRO ECHAV~ BUENOS AIRES 1094 0.00640
## 5 ALFONSO LÓPEZ CASTILLA 1994 0.00752
## 6 ALTAMIRA ROBLEDO 755 0.00927
## 7 ALTAVISTA ALTAVISTA 168 0.0119
## 8 ALTAVISTA BELEN 1988 0.00553
## 9 ALTAVISTA SECTO~ ALTAVISTA 507 0.00592
## 10 ALTOS DEL POBLA~ EL POBLADO 375 0.00267
## # ... with 298 more rows, and 6 more variables: proporcion_pension <dbl>,
## # subsidio <dbl>, salario <dbl>, arriendo <dbl>, total_gastos <dbl>,
## # gasto_servicios <dbl>Escritura.
write_excel_csv2(db_ingresos, "databases/db_ingresos.csv")
write_excel_csv2(db_ingresos[complete.cases(db_ingresos),], "databases/db_ingresos_completecases.csv")Esta dimensión llamada “ingresos” aborda algunas preguntas referentes a los ingresos y gastos de las personas de los diferentes barrios de Medellín.
DESCRIPCIÓN DE LAS VARIABLES:
- encuesta_calidad_barrio:Barrios de Medellín.
- encuesta_calidad_comuna:Comunas de Medellín.
- n:Cantidad de personas encuestadas por cada barrio.
- proporcion_ingreso_extra:Proporción de personas por cada barrio que reciben ingreso sin tener que trabajar por parte de algún trabajo o negocio en los diferentes barrios.
- proporcion_pension:Proporción de personas por cada barrio que reciben pensión o jubilación.
- subsidio:Promedio por cada barrio de cuánto recibió en promedio mensual la persona durante los últimos 12 meses por concepto de ayudas en dinero, intereses o dividendos u otras fuentes como auxilio o subsidio de transporte, subsidio familiar y subsidio educativo.
- salario:Promedio por cada barrio de cuánto recibió la persona el mes pasado de haber hecho la encuesta por concepto de trabajo.
- arriendo:Pago promedio en cada barrio por: Pago(Vivienda en arriendo o subarriendo mensual) o Cuota (Propia, la están pagando mensual).
- total_gastos:Promedio por cada barrio de cuánto asciende el total de gastos mensuales en cada hogar.
gasto_servicio:Es el promedio por cada barrio de la cantidad del ingreso mensual destina a servicios públicos en cada hogar.
Nota:En la base de datos utilizada se encontraron algunas celdas vacías en las preguntas utilizadas para calcular las variables “subsidio”, “salario”, y “arriendo”, esto puede deberse a que la persona encuestada no recibe algún tipo de subsidio o salario, y no paga arriendo. Al fijarnos en los barrios en los que para ninguna persona hay respuesta, nos dimos cuenta que son barrios de invasión, por lo cual es coherente que haya un cero para tal respuesta; también se pudo observar que en algunos casos las personas no reciben sueldo, pero sí reciben un subsidio. Por ésto, en las celdas de nuestra nueva base de datos db_ingresos en donde no fue posible calcular promedios, pues el cociente quedaba en cero y esto tiende a infinito, se reemplazó por ceros.
Análisis descriptivo
summary(db_ingresos)## encuesta_calidad.barrio encuesta_calidad.comuna n
## Length:308 Length:308 Min. : 4.0
## Class :character Class :character 1st Qu.: 348.8
## Mode :character Mode :character Median : 862.5
## Mean :1073.2
## 3rd Qu.:1439.0
## Max. :8987.0
## proporcion_ingreso_extra proporcion_pension subsidio
## Min. :0.000000 Min. :0.1000 Min. : 0
## 1st Qu.:0.003061 1st Qu.:0.4049 1st Qu.: 244619
## Median :0.005511 Median :0.4426 Median : 370100
## Mean :0.006463 Mean :0.4386 Mean : 520659
## 3rd Qu.:0.008005 3rd Qu.:0.4696 3rd Qu.: 587997
## Max. :0.066667 Max. :0.7500 Max. :8241967
## salario arriendo total_gastos gasto_servicios
## Min. : 0 Min. : 0 Min. : 425763 Min. : 50500
## 1st Qu.: 305632 1st Qu.: 283998 1st Qu.: 833227 1st Qu.:150364
## Median : 437404 Median : 409747 Median :1087169 Median :206197
## Mean : 703446 Mean : 535185 Mean :1393247 Mean :225134
## 3rd Qu.: 742865 3rd Qu.: 654791 3rd Qu.:1580345 3rd Qu.:275162
## Max. :8734250 Max. :1901874 Max. :4976319 Max. :588615- Dada la información anterior, se observa que en general son pocas las personas que tienen algún tipo de ingreso por algún trabajo o negocio sin tener que trabajar, habiendo barrios en los que la proporción de estas personas es cero, es decir, ninguna persona recibe ingresos sin tener que trabajar. Boquerón es el barrio donde la proporción es mayor. Por otro lado vemos que en todos los barrios hay personas que reciben pensión o jubilación, siendo Guayaquil en donda está la mayoría de personas y San José en donde está la minoría.
- Se logra ver que hay barrios en los que ninguna persona recibió durante los últimos 12 meses concepto de ayudas en dinero, intereses o dividendos u otras fuentes como auxilio o subsidio de transporte, subsidio familiar y subsidio educativo y el barrio.
Como se hablaba anteriormente, hay barrios en los que las personas para el tiempo en que fue hecha la encuesta no recibieron un salario el mes anterior, por el contrario en los barrios en donde sí reciben salarios, Villa de Aburrá es el que tiene un promedio más alto (8.734.250).
- Haciendo una comparación del pago de arriendo promedio, el barrio de Las Lomas No.2 es donde se realiza un mayor pago y el barrio en donde hay un mayor promedio de gastos totales por hogar es en el Tesoro y cerca de ahí por los Altos del Poblado es el barrio en donde hay un mayor promedio de dinero destinado para el pago de servicios públicos.
min_gasto_servicios<-apply(db_ingresos[10],2,min, na.rm=T)
db_ingresos[which(db_ingresos$gasto_servicios==min_gasto_servicios),1]## # A tibble: 1 x 1
## encuesta_calidad.barrio
## <chr>
## 1 SAN JOSÉmax_gasto_servicios<-apply(db_ingresos[10],2,max, na.rm=T)
db_ingresos[which(db_ingresos$gasto_servicios==max_gasto_servicios),1]## # A tibble: 1 x 1
## encuesta_calidad.barrio
## <chr>
## 1 ALTOS DEL POBLADOGráfico (1).
grafico<-ggpairs(db_ingresos,columns = 4:10, mapping=aes(alpha = 0.4))
ggsave(filename="images/correlacion.png", plot =grafico, width = 20, height = 20)
Correlación entre las variables
Matriz de correlaciones:
data<-as.matrix(db_ingresos[,3:9])
mcor<-cor(data)
mcor## n proporcion_ingreso_extra
## n 1.00000000 -0.06518933
## proporcion_ingreso_extra -0.06518933 1.00000000
## proporcion_pension -0.17011586 -0.03081172
## subsidio -0.01152908 -0.02913350
## salario -0.02949656 -0.01847511
## arriendo -0.10414827 0.01745026
## total_gastos -0.10706213 0.01250525
## proporcion_pension subsidio salario
## n -0.17011586 -0.01152908 -0.02949656
## proporcion_ingreso_extra -0.03081172 -0.02913350 -0.01847511
## proporcion_pension 1.00000000 0.03797035 0.08478863
## subsidio 0.03797035 1.00000000 0.28460736
## salario 0.08478863 0.28460736 1.00000000
## arriendo 0.05240449 0.40708611 0.65315248
## total_gastos 0.02899374 0.40252794 0.63576569
## arriendo total_gastos
## n -0.10414827 -0.10706213
## proporcion_ingreso_extra 0.01745026 0.01250525
## proporcion_pension 0.05240449 0.02899374
## subsidio 0.40708611 0.40252794
## salario 0.65315248 0.63576569
## arriendo 1.00000000 0.95678984
## total_gastos 0.95678984 1.00000000Gráfico (2) Correlaciones:
corrplot(mcor,type = "upper",order = "hclust",tl.col = "black",tl.srt = 45)
- En los gráficos (1) y (2) se observa que hay un par de variables que presentan una alta correlación, que son las correspondietes al pago de arriendo y al total de gastos por hogar. Otras variables que presentan una correlación leve son las referentes al salario recibido con el pago de arriendo y al salario recibido con el total de gastos; estas correlaciones eran de esperarse.
Escalemos.
db_ingresos_escalado <- db_ingresos%>%
mutate_at(vars(-n, -encuesta_calidad.barrio, -encuesta_calidad.comuna), scale)Selección de las columnas para el agrupamiento.
datos_escalados_agrupamiento <- db_ingresos_escalado %>%
dplyr::select(-n, -encuesta_calidad.barrio, -encuesta_calidad.comuna)Apliquemos análisis de componente principales (PCA).
pca <- prcomp(datos_escalados_agrupamiento)
summary(pca)## Importance of components:
## PC1 PC2 PC3 PC4 PC5 PC6 PC7
## Standard deviation 1.8902 1.0168 0.9845 0.8853 0.68540 0.3579 0.20592
## Proportion of Variance 0.5104 0.1477 0.1385 0.1120 0.06711 0.0183 0.00606
## Cumulative Proportion 0.5104 0.6581 0.7966 0.9085 0.97564 0.9939 1.00000Varianza explicada por cada componente
porcentaje_exp <- 100 * (pca$sdev ^ 2) / sum(pca$sdev ^ 2)
ggplot(data = data.frame(componente = factor(1:ncol(datos_escalados_agrupamiento)), porcentaje_exp)) +
geom_col(mapping = aes(x = componente, y = porcentaje_exp))
suma_acumulada_porcentaje <- cumsum(porcentaje_exp)
ggplot(data = data.frame(componente = factor(1:ncol(datos_escalados_agrupamiento)), suma_acumulada_porcentaje),
aes(x = componente, y = suma_acumulada_porcentaje, group = 1)) +
geom_point() +
geom_line() +
geom_hline(yintercept = 90, color = "red") +
geom_vline(xintercept = 7)
suma_acumulada_porcentaje## [1] 51.03878 65.80981 79.65645 90.85302 97.56407 99.39425 100.00000- Notemos que con las primeras cuatro componente se explica la variablidad en un 90.8 %
Proyección respecto a las componentes principales
datos_proyectados <- t(t(pca$rotation) %*% t(datos_escalados_agrupamiento))datos_agrupamiento <- datos_proyectados[, 1:4] %>% as_tibble()Agrupamiento mediante K-means
Se utilizará la metodología de Elbow para identificar el k óptimo que ofreza una menor distancia entre clusters como medida de error Wss.
datos_escalados_agrupamiento<- db_ingresos_escalado %>%
dplyr::select(-n, -encuesta_calidad.barrio, -encuesta_calidad.comuna)
k_maximo <- 15
wss <- sapply(1:k_maximo,
function(k) {
kmeans(datos_escalados_agrupamiento, k, nstart = 10, iter.max = 15 )$tot.withinss
})
plot(1:k_maximo,wss,
type="b", pch = 19, frame = FALSE,
xlab="Número de grupos [K]",
ylab="Total de suma de cuadrados intra-grupos")
Según nuestro algoritmo de K-means, nuestro K óptimo es de 9 grupos, pero se considera más pertinente ultilizar un número menor para favorecer la interpretabilidad.
k_optimo_kmeans<-9
k_seleccionado<-5Agrupamiento con K óptimo
Fijaremos una semilla de los número aleatorias para garantizar la reproducibilidad de las conclusiones del agrupamiento.
set.seed(3)
agrupamiento_kmeans <- kmeans(datos_agrupamiento, centers = k_seleccionado, nstart = 10, iter.max = 15 )
db_ingresos_kmeans <- bind_cols(db_ingresos, cluster = agrupamiento_kmeans$cluster)Mapa geográfico mediante el agrupamiento K-means
political <- shapefile("Barrio_Vereda/Barrio_Vereda.shp")
Encoding(political@data$NOMBRE) <- "UTF-8"
political$NOMBRE <- political$NOMBRE %>% toupper() %>% str_replace("DE MESA", "DE MESA")
grupos_barrios <- data.frame(barrio_nombre = db_ingresos$encuesta_calidad.barrio, grupo = agrupamiento_kmeans$cluster)
nombres_mapa <- data.frame(nombre_barrio = political$NOMBRE)
vector_nombres = c()
vector_grupos = c()
for(nombre_mapa in nombres_mapa$nombre_barrio) {
grupo <- grupos_barrios[grupos_barrios$barrio_nombre == nombre_mapa, 2][1]
vector_nombres <- c(vector_nombres, nombre_mapa)
vector_grupos <- c(vector_grupos, grupo)
}
factpal <- colorFactor(rainbow(k_seleccionado), vector_grupos)
mapa_grupos <-tibble(
nombre_barrio = vector_nombres,
grupo = vector_grupos,
color = factpal(grupo),
descripcion = case_when(
grupo == 1 ~ "Grupo 1",
grupo == 2 ~ "Grupo 2",
grupo == 3 ~ "Grupo 3",
grupo == 4 ~ "Grupo 4",
grupo == 5 ~ "Grupo 5",
TRUE ~ "Sin muestra representativa")
)
colores <- mapa_grupos %>%
dplyr::select(-nombre_barrio) %>%
distinct(color, .keep_all = TRUE) %>%
arrange(desc(grupo))
leaflet(data = political) %>%
addTiles() %>%
addProviderTiles(providers$OpenStreetMap) %>%
addPolygons(fill = TRUE, stroke = TRUE, weight = 2, color = mapa_grupos$color,
label = as.character(political$NOMBRE),
popup = as.character(political$NOMBRE)) %>%
addLegend("bottomright", colors = colores$color, labels = as.character(colores$descripcion))